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Abstract 

Novel  parametric  reduced-order  models  are  proposed  for  fast  reanalysis  to  pre¬ 
dict  the  dynamic  response  of  complex  structures,  which  suffered  thickness  vari¬ 
ations  caused  by  design  changes  or  damage  in  one  or  more  substructures.  Para¬ 
metric  reduced-order  models  developed  previously  have  two  important  challenges 
to  overcome  to  improve  accuracy  and  performance:  (a)  the  transformation  matrix 
is  not  mathematically  stable,  (b)  the  Taylor  series  parameterization  techniques  do 
not  capture  thickness  variations  of  the  structure  modeled  with  solid-type  elements 
due  to  the  highly  nonlinear  dependence  on  thickness  changes.  Thus,  herein,  a 
new  transformation  matrix  and  novel  parameterization  techniques  are  proposed. 
Usual  reduced-order  models  have  an  additional  challenge,  namely  the  difficulty  in 
reducing  the  interface  degrees  of  freedom.  Thus  a  way  of  reducing  the  interface 
degrees  of  freedom  is  also  proposed.  The  predicted  vibration  responses  of  com¬ 
plex  structures  are  shown  to  agree  very  well  with  results  obtained  using  a  much 
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1.  Introduction 

Component  mode  synthesis  techniques  are  well  established  [1,  2,  3,  4,  5,  6,  7] 
in  the  field  of  structural  dynamic  analysis  as  an  alternative  to  conventional  finite 
element  models  (FEMs)  with  large  number  of  degrees  of  freedom  (DOFs).  Com¬ 
ponent  mode  synthesis  (CMS)  belongs  to  a  wide  class  of  domain  decomposition 
techniques.  CMS  is  a  substructural  based  technique,  which  divides  the  global 
structure  into  several  substructures,  and  the  DOFs  of  those  substructures  are  re¬ 
duced  significantly.  Then,  each  individual  substructure  in  the  CMS  domain  is 
reconnected,  and  the  system  dynamic  responses  are  predicted  very  efficiently  and 
accurately. 

CMS  has  become  a  very  popular  numerical  tool  in  aerospace  and  automotive 
engineering  because  it  usually  meets  high  standards  of  computational  efficiency. 
Computational  efficiency  is  illustrated  by  significant  cost  saving  when  remeshing 
is  needed,  since  this  task  can  be  done  locally,  i.e.  on  each  substructure  separately. 
However,  the  remeshing  process  might  also  be  time  consuming  computationally 
and  manually  for  design  purposes  such  as  structural  optimization,  and  for  damage 
modeling  for  structural  health  monitoring.  Therefore,  reduced-order  modeling 
techniques  for  design  and  damage  modeling  purposes  are  needed. 

The  reduced-order  models  (ROMs)  for  design  and  damage  modeling  were  in¬ 
troduced  almost  fifteen  years  ago  by  Balmes  et  al.  [8,  9]  to  avoid  the  relatively 
expensive  process  of  reanalysis  of  complex  structures.  In  addition,  several  other 
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ROMs  referred  to  as  parametric  reduced-order  models  (PROMs)  have  been  devel¬ 
oped  [10,  11, 12,  13].  In  particular,  multi-component  PROMs  (MC-PROMs)  have 
been  developed  recently  by  Hong  et  al.  [13].  For  robust  substructure  (re)analysis, 
MC-PROMs  are  advantageous  because  they  allow  several  substructures  to  have 
parametric  variability  in  characteristics  such  as  geometric  parameters  (e.g.,  thick¬ 
ness),  or  material  properties  (e.g.,  Young’s  modulus).  MC-PROMs  are  perfectly 
suited  for  predicting  the  vibration  response  of  structures  modeled  with  shell-type 
finite  elements  which  can  have  thickness  variations.  However,  if  the  structure  is 
modeled  with  brick-type  finite  elements,  and  if  the  brick-type  elements  require 
local  volume  changes  during  reanalysis,  then  MC-PROMs  cannot  be  effectively 
used  to  predict  the  dynamic  response.  This  is  because  MC-PROMs  use  third-order 
Taylor  series  for  parameterization.  These  Taylor  series  do  not  capture  accurately 
the  variation  of  the  mass  and  stiffness  matrices  for  brick-type  finite  elements  be¬ 
cause  the  volume  of  local  finite  elements  can  change  during  the  reanalysis.  Con¬ 
sequently,  some  entries  of  the  mass  and  stiffness  matrices  for  brick-type  finite 
elements  vary  highly  nonlinearly  with  respect  to  geometric  variations  in  the  struc¬ 
ture.  Herein,  a  novel  parameterization  technique  is  proposed  to  capture  these 
element-level  nonlinearities. 

Another  challenge  for  MC-PROMs  is  that  they  can  be  numerically  not  stable 
due  to  the  transformation  matrix  they  employ.  Specifically,  the  transformation 
matrix  consists  of  static  constraint  modes  and  fixed  interface  normal  modes  com¬ 
puted  for  a  set  of  nominal  parameters,  and  a  few  sets  of  perturbed  parameter  values 
(typically  up  to  3  sets  per  parameter)  [13].  If  all  static  constraint  modes  are  kept 
and  many  normal  modes  are  included,  then  the  size  of  the  system-level  mass  and 
stiffness  matrices  can  be  nearly  singular  (and  can  even  be  larger  than  that  of  the 
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full-order  models).  This  is  because  the  transformation  matrix  can  contain  vectors 
which  are  nearly  linearly  dependent.  These  vectors  are  usually  normal  modes  for 
substructures  where  the  parametric  variation  (e.g.,  modulus  of  elasticity)  does  not 
affect  the  component-level  normal  modes.  Moreover,  the  transformation  matrix 
used  in  MC-PROMs  was  designed  for  small  parameter  variations  which  ensures 
that  the  space  spanned  by  the  basis  vectors  at  a  component-level  does  not  depend 
nonlinearly  on  parameter  variations.  That  approach  can  break  down  because  of 
the  volume  variations  which  can  occur  when  brick  elements  are  used. 

Another  challenge  of  CMS  methods  and  MC-PROMs  is  that  they  often  re¬ 
quire  an  excessively  large  number  of  interface  DOFs  because  (often)  these  DOFs 
are  many  and  are  hard  (or  impossible)  to  reduce.  To  address  this  issue,  Castanier 
et  al.  [6]  proposed  that  the  physical  interface  DOFs  be  replaced  by  global  interface 
modes,  which  were  also  called  characteristic  constraint  (CC)  modes.  However, 
this  concept  is  not  optimal  for  substructural-based  techniques  because  CC  modes 
are  system-level  interface  modes,  not  sub  structural- level  interface  modes.  Thus, 
a  new  technique  to  reduce  the  interface  DOFs  locally  is  proposed  herein  and  re¬ 
ferred  to  as  local-interface  reduction. 

This  paper  is  organized  as  follows.  In  Section  2,  the  element-level  nonlinear¬ 
ity  due  to  the  volume  variations  of  finite  elements  of  brick  or  other  type  is  evalu¬ 
ated,  and  a  novel  parameterization  technique  is  proposed  to  capture  this  nonlinear¬ 
ity.  Next,  in  Section  3,  CMS  is  briefly  reviewed,  and  next-generation  parametric 
reduced-order  models  (NX-PROMs)  are  proposed.  In  Section  4,  to  locally  reduce 
the  interface  DOFs,  a  local-interface  reduction  technique  is  presented.  Section 
5  discusses  the  procedure  to  assemble  substructural  mass  and  stiffness  matrices 
(with  and  without  implementing  the  local-interface  reduction  technique).  In  Sec- 
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tion  6,  numerical  examples  such  as  a  plate  structure,  an  L-shaped  structure,  and  a 
realistic  vehicle  model  (a  high  mobility  multipurpose  wheeled  vehicle,  HMMWV) 
modeled  with  brick-type  finite  elements  are  used  to  demonstrate  the  proposed 
methods.  Finally,  conclusions  are  summarized  in  Section  7. 

2.  Robust  Parameterization  Techniques  for  Element-Level  Nonlinearity 

For  structural  design  and  damage  modeling  purposes,  the  parameterization  of 
the  mass  and  stiffness  matrices  can  be  the  most  important  step.  This  is  because  the 
parameterization  techniques  enable  capturing  mass  and  stiffness  variations  due  to 
design  changes  or  damage  in  the  structure.  Thus,  the  reanalysis  time  can  be  sig¬ 
nificantly  reduced  because  the  finite  element  mesh  does  not  need  to  be  modified 
and  remodeled.  The  parameterization  technique  has  to  be  adapted  for  the  differ¬ 
ent  characteristics  of  each  type  of  finite  element  used.  For  example,  the  thickness, 
Young’s  modulus,  and  material  density  variations  of  shell-type  elements  can  be 
captured  well  by  third-order  Taylor  series  [8,  9,  12,  13].  However,  we  found 
that  the  thickness  variations  for  a  brick  and  other  types  of  finite  elements  such 
as  hexagonal  and  tetrahedron  elements  cannot  be  captured  well  by  Taylor  series 
of  low  order  due  to  an  element-level  nonlinear  characteristics  caused  by  volume 
variations.  For  elements  which  have  volume,  local  thickness  variations  induce 
volume  variations  in  the  elements.  In  contrast,  Taylor  series  works  well  for  pa¬ 
rameterizing  shell-type  elements  because  these  do  not  have  actual  volume.  In  this 
section,  a  parameterization  technique  that  captures  thickness  variations  of  brick 
and  other  types  of  finite  elements  is  the  focus. 

Fig.  1  shows  an  8-node  brick-type  element  which  uses  first-order  (linear)  shape 
functions.  Coordinates  x,  y,  and  z  are  global,  and  coordinates  £,  r/,  and  (  are  local. 
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As  a  conceptual  example,  consider  that  the  four  nodes  on  the  top  surface  in  Fig.  1 
move  by  a  distance  At.  The  brick-type  element  has  a  volume,  so  when  each  node 
moves,  the  volume  of  the  brick-type  element  varies.  Thus,  the  parameterization 
technique  has  to  account  for  these  volume  variations.  To  that  aim,  let  us  first  revisit 
the  formulation  used  to  derive  stiffness  matrices  for  brick-type  elements  [14,  15]. 
The  equation  used  to  obtain  the  stiffness  matrix  can  be  expressed  as 


K 


B  BBdV 


' v 


<=i  rv= 1  A=1 


BTBBd{dr}d( 


'<=— 1  J  )7=-l  «/£=-! 
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EES  WiWjWkBT(£i,  rjj,  (fc)DB((j,  ty,  Cfc)det(J (&,  ry,  ( '*)), 

i= 1  j= 1  k= 1 


where  B  is  a  strain  matrix  (which  contains  derivatives  of  the  linear  shape  functions 
in  global  coordinates),  and  D  is  an  elasticity  matrix  (which  contains  Poisson’s  ra¬ 
tio  v,  and  the  elastic  modulus  E ).  The  determinant  of  the  Jacobian  in  Eq.  (1)  is 
obtained  from  the  coordinate  transformation  of  the  strain  matrix  B.  The  determi¬ 
nant  contains  in  its  denominator  a  cubic  polynomial  of  £,  rj  and  (,  which  reflects 
volume  variations.  Thus,  the  parameterization  should  also  contain  a  cubic  poly¬ 
nomial  in  the  denominator.  To  establish  the  coefficients  of  this  cubic  polynomial, 
the  volume  variations  of  brick-type  elements  are  considered.  As  shown  in  Fig.  1, 
one  or  several  nodes  on  the  top  surface  move  to  capture  thickness  variations.  The 
resulting  volume  can  be  expressed  as 


V  =  Vn 


1  +  d 


P~P  o 
Po 


—  Vo  1  +  d 


Ap\ 


Po  ) 


where  V  and  Vo  are  the  final  and  the  initial  volume  of  the  brick-type  element, 
and  p  and  p0  are  the  target  parameter  value  (thickness)  and  the  initial  parameter 
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value.  When  only  one  node  on  the  top  surface  moves,  the  coefficient  d  is  1/3. 
Similarly,  when  two  nodes  move,  d  =  1/2.  Also,  when  three  nodes  move,  d  —  1. 
Finally,  when  four  nodes  move,  the  volume  variation  is  proportional  to  A p,  i.e. 
V  =  I  ;,  =^- .  This  last  type  of  variation,  proportional  to  A p,  is  very  well  captured 
by  a  regular  interpolation.  The  other  three  cases,  however,  are  not.  To  address 
this  issue,  a  cubic  polynomial  which  considers  the  volume  variation  of  brick-type 
elements  with  a  target  parameter  variation  A p  is  defined  as 


D(Ap)  =  (  1  +  Af 
Po  J 


1 A  p\ 


1  A  p\ 


1  +  -—  1  (1  +  -—  • 

3  Po  ) 
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2  Po  / 

The  new  parameterization  equation  consists  of  a  fourth-order  interpolation  in  the 
numerator  and  the  cubic  polynomial  in  Eq.  (2)  in  the  denominator,  which  yields 
the  new  parameterization  as 

K0  +  KXA  p  +  K2(Ap)2  +  K3(Ap)3  +  K4(Ap)4 


K(p0  +  A p) 


D(Ap) 


(3) 


To  calculate  the  matrices  K0,  K4,  K2,  K3  and  K4  in  Eq.  (3),  five  equations  are 
needed.  For  that,  stiffness  matrices  for  five  parameter  values  are  computed.  First, 
consider  the  case  where  A p  =  0.  One  obtains 


K(p0)  «  K0. 


(4) 


Next,  consider  A p  =  iSp  (with  i  =  1,  2,  3, 4),  one  obtains 


K(p0  +  iSp)  « 


K0  +  K  ^iSp)  +  K2(i5p)2  +  K  3(iSp)3  +  K  4(iSp)4 
D(Sp) 


(5) 


Rearranging  Eqs.  (4)  and  (5)  into  matrix  form,  for  each  entry  e,q  of  the  matrices 
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Ki,  K2,  K3  and  K4,  one  obtains 
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Equation  (6)  can  be  easily  solved  by  simply  inverting  the  5x5  matrix  on  the  left 
hand  side.  This  matrix  is  non-singular  and  very  well  behaved  for  inversion.  Also, 
note  that  this  inversion  has  to  be  done  only  once  (for  a  given  Sp).  Let  us  denote 
by  A  this  inverse  matrix,  i.e. 
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Re-arranging  Eq.  (6)  using  the  entries  in  A,  one  obtains 


K(p0  +  Ap)  ft*  (An  +  A21Ap  +  A31Ap2  +  A41Ap3  +  A51Ap4)K(p0) 

+  (M2  +  A22Ap  +  A32Ap2  +  A42Ap3  +  M2Ap4)K(i°o  +  Sp) 

+  (-4 13  +  M3  A  p  +  A33Ap2  +  A43Ap3  +  A53Ap4)K(p0  +  2  Sp) 

+  (A14  +  A24Ap  +  A34Ap2  +  A44Ap3  +  A34Ap^)lA{pQ  +  3  Sp) 

+  (Ms  +  MsAp  +  A35Ap2  +  A45Ap3  +  MsAp4)K(p0  +  4 Sp),  or 


K(p0  +  A p)  « 


50K(p0)  +  &iK(p0  +  Sp)  +  &2K(p0  +  25p) 


+  53K(p0  +  3Sp)  +  b4K(p0  +  Up). 


(8) 
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Equation  (8)  shows  that  K(p0  +  A p)  is  simply  a  linear  combination  of  five  (pre¬ 
computed)  matrices.  The  coefficients  in  the  linear  combination  depend  very  non- 
linearly  on  A p.  That  is  the  key  factor  which  ensures  the  high  accuracy  of  the  new 
parameterization.  Note  that  the  computational  cost  of  the  new  parameterization  is 
the  same  as  that  for  a  regular  fourth-order  interpolation.  The  accuracy,  however, 
is  higher  (as  shown  on  Section  6.1). 

3.  Reduced-Order  Models 

The  parameterization  techniques  proposed  in  Section  2  are  for  the  full-order 
finite  element  model.  However,  the  main  objective  of  this  work  is  to  predict  vi¬ 
bration  responses  using  ROMs  (as  opposed  to  full-order  models)  to  reduce  the 
calculation  time.  To  detail  the  construction  of  ROMs,  the  fixed-interface  Craig- 
Bampton  component  mode  synthesis  (CB-CMS)  [5]  is  reviewed  briefly.  Next, 
a  new  transformation  matrix  is  presented  and  used  in  conjunction  with  the  new 
parameterization  technique  discussed  in  Section  2.  Finally,  NX-PROMs  are  con¬ 
structed. 

3.1.  Brief  review  of  Craig -Bampton  component  mode  synthesis 

In  this  section,  the  fixed-interface  CB-CMS  [2]  method  is  reviewed.  This 
modeling  approach  is  broadly  used  because  of  its  simplicity  and  computational 
stability.  To  apply  the  CB-CMS,  the  complex  structure  of  interest  is  partitioned 
into  substructures.  The  DOFs  of  each  substructure  are  further  partitioned  into 
active  DOFs  on  the  interface  (indicated  by  the  superscript  A),  and  omitted  DOFs 
in  the  interior  (indicated  by  the  superscript  O).  The  mass  and  stiffness  matrices 
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for  a  component  i  can  then  be  partitioned  to  obtain 


mAA  mAO 

’kAA 

1 

o 

L  l 

,  and  Kj  = 

l 

m,fA  mf° 

_k?A 

koo  j 

Next,  the  physical  coordinates  are  changed  to  a  set  of  coordinates  representing 
the  amplitudes  of  a  selected  set  of  fixed-interface  component-level  normal  modes 
(I» ,v  (indicated  by  the  superscript  N ),  and  the  amplitudes  of  the  full  set  of  static 
constraint  modes  =  — kf°  1kfA  (indicated  by  the  superscript  C).  The  trans¬ 
formed  mass  and  stiffness  matrices  for  component  i  can  be  expressed  as 


y^CBCMS 


~  r1 

my 


m.?N 


m 


NC 


m 


N 


and 


j^CBCMS 


k? 

kP 

kfc 

kf 

In  this  work,  the  CB-CMS  method  is  used  only  for  the  substructures  which  do 
not  have  any  parameter  variation  or  damage. 


3.2.  Next- generation  parametric  reduced-order  models 

In  this  section,  MC-PROMs  are  improved  to  be  more  robust  and  mathemati¬ 
cally  stable.  The  resulting  models  are  referred  to  as  NX-PROMs. 


3.2.1.  Transformation  matrix  for  NX-PROMs 

The  transformation  matrix  for  NX-PROMs  is  constructed  somewhat  similar 
to  MC-PROMs  (which  was  constructed  by  using  the  idea  behind  CB-CMS).  It 
also  has  a  set  of  static  constraint  modes  and  a  set  of  fixed-interface  normal 
modes  d> A .  However,  the  transformation  matrix  for  NX-PROMs  has  a  different 
set  of  static  constraint  modes  and  a  different  set  of  fixed-interface  normal  modes 
compared  at  CB-CMS  and  MC-PROM.  This  transformation  matrix  can  be  written 
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for  component  i  as 


I  0 

q>C  $>N 

aug,2  aug,2 


where  is  referred  to  as  the  matrix  of  augmented  constraint  modes 


\T 

aug,i 


^  0,2 


c 


^  4,2 


5 


and  is  referred  to  as  the  matrix  of  augmented  fixed-interface  normal  modes 


aug,2 


N 


H 


H 


4,2 


Matrices  \I7p7  and  <1>(Y  correspond  to  the  nominal  parameter  values,  whereas  ma¬ 
trices  SPp  and  d>;Y  (7  =  1,  2,  3, 4)  correspond  to  four  other  parameter  values. 

In  general,  the  columns  of  are  not  orthogonal.  Therefore,  for  numerical 
stability,  an  orthogonal  basis  for  the  space  spanned  by  these  modes  is  computed. 
This  basis  is  obtained  by  a  truncated  set  of  left  singular  vectors  U;Y  of  [13]. 
Thus,  the  new  transformation  matrix  can  be  expressed  as 


T,  = 


c 

aug,2 


o 

uf 


(9) 


The  transformation  matrix  in  Eq.  (9)  can  be  used  to  project  the  physical  domain 
onto  the  NX-PROM  domain.  The  stiffness  matrix  K;YXPROM  =  T’/  Kf  p  ,,/T,  for 
component  i  (where  Kf  EM  represents  the  stiffness  matrix  of  the  full-order  model 
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of  component  i)  can  be  partitioned  to  obtain 


K 


NXPROM 


K 

K 

K 

K 

K 
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c 

00,  i 

Kc 

■^oi,* 

K c 

^02,1 
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■*^03,* 

Kc 
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kcn 
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C 
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Kc 

-■-Ml,* 

Kc 

^12,1 

Kc 

-*^13  ,i 

Kc 

-*^14,* 

Kcn 

^ii,* 

C 

20,  * 

Kc 

^21,  i 

K  c 

*^22,1 

Kc 

■rv23,* 

Kc 

J:v24,* 

Kcn 

*^22  ,i 

C 

30,  * 

Kc 

*^31,1 

Kc 

iv32,  i 

KC 

*^33, i 

Kc 

^34,1 

kcn 

■*^33,* 

C 

40,* 

Kc 

Kc 

^42,* 

KC 

Kc 

J:v44,* 

NC 

00,* 

knc 

knc 

*^22, i 

knc 

*^33  ,i 

knc 

J:v44,* 

Kf 

A  similar  relation  is  obtained  for  the  mass  matrix  of  component  i.  Here,  the 
DOFs  corresponding  to  the  constraint  part  (superscript  C )  are  repeated  for  the 
five  parameter  values  (denoted  by  subscript  0,  1,  2,  3  and  4).  Note  that  in  such  an 
approach,  the  size  of  the  mass  and  stiffness  matrices  can  be  quite  large.  Also,  these 
matrices  may  be  ill-conditioned  because  the  columns  of  ^  are  not  necessarily 
linearly  independent. 

To  address  this  issue,  a  new  method  to  account  for  the  static  constraint  modes 
is  developed.  This  new  method  avoids  duplicating  the  interface  DOFs  (C)  and 
captures  the  interface  effects  more  accurately.  The  approach  reduces  the  number 
of  sets  of  static  constraint  modes  (used  to  obtain  MfrxpROM  and  Kf XPROM) 
from  five  sets  to  just  one  set.  Consider  that  the  actual  value  p  of  the  parameter 
where  the  reanalysis  is  needed  exists  between  the  Ith  and  (7  +  l)th  parameter 
values  (/  =  0, 1,  2,  3  or  4)  which  were  used  to  construct  ~MfXPROA1 .  In  that  case 
(described  in  Fig.  2),  a  new  static  constraint  mode  can  be  generated  by  linearly 
interpolating  between  the  static  constraint  modes  for  the  Ith  and  the  (/  +  l)th 
parameter  values  to  obtain 


f  Pi+ i  ~p\ 
\Pi+i-Pi) 


f  P~Pi  \ 
\Pi+ 1  ~  Pi) 


\l> 


c  _ 

l+l,i  ~ 


c 

l,i 


C 


(10) 
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This  new  static  constraint  mode  T'f  replaces  all  five  static  constraint  modes  used 
to  construct  T,  in  Eq.  (9). 

Note  that  this  reduction  can  be  implemented  without  re-constructing  T\:  for 
each  case  of  parameter  variation.  Instead,  only  simple  linear  combinations  of 
partitions  of  the  matrices  m nxprom  anq  r nxprom ,  are  nee(je(j  Details  are 
given  in  Section  3.2.2  below.  In  the  end,  the  final  NX-PROM  mass  and  stiffness 
matrices  have  only  a  single  set  of  constraint  modes  'Pf  which  always  has  linearly 
independent  columns. 


3.2.2.  parameterization  for  NX-PROMs 

The  new  interpolation  presented  in  Section  2  is  applied  to  NX-PROMs.  To 
that  aim,  five  mass  and  five  stiffness  are  constructed  for  each  component  i  (for 
l  —  0,1, 2, 3, 4)  as  follows 

Mxx  =  TjM(p0  +  Up) fi  and  Kxx  =  fjK(p0  +  Up) TV 


These  matrices  are  not  all  used  independently  to  form  NX-PROMs.  Instead,  they 
are  linearly  combined  to  implement  the  single  set  of  static  constraint  modes  vPp 
in  Eq.  (10).  Thus,  conceptually,  T*  is  replaced  by  TV  given  by 


I  0 


q>c  ttat 
aug,z  i 


(cTjR^j  T  /^iRj+l^) 


T  (cTjR,;  j  T  /5jR/-|-l,j)  j 


where  and  R./  h ,  a  are  masking  matrices  of  zeros  and  ones.  The  matrix  R,;.,  is 

a  6  x  2  block  matrix  where  the  blocks  (/,  1)  and  (6,  2)  are  unit  matrices,  while  all 
other  blocks  are  zero.  The  first  five  rows  correspond  to  (I  =  0, 1,  2,  3, 4)  and 
the  last  row  corresponds  to  U^.  This  new  transformation  matrix  T,  is  applied  to 
the  mass  and  stiffness  matrices  of  component  i  to  construct  NX-PROMs.  First, 
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five  mass  and  stiffness  matrices  are  constructed  for  each  parameter  variation  s5p 
(for  s  —  0, 1,  2,  3, 4)  as  follows 


M™  =  (q,RT  +  p, RT  M)  M™  (o,R,f  +  ftRi+1,i) , 

K™  =  (<*,R£  +  AR£m)  Kj?  (<*R«  +  m*u)  ■ 

Next,  the  new  interpolation  discussed  in  Section  2  is  applied  using  and 


Eq.  (8)  is  used  to  obtain 


S,l 


M(p0  +  A  p){ 


NX 


bo  M™  +  b,  Miy  +  b2mlt  +  63MAA  +  64m^ 


\NX 


NX 


tNX 


iNX 


K(p0  +  A p) 


NX 


b0  +  6iK  AA  +  ^KAA  +  63KAA  +  64K 


-jvx 


NX 


-NX 


-NX 

L4,i 


where  bs  (s  =  0, 1,  2,  3, 4)  are  computed  for  each  A p  by  using  the  matrix  A  and 
the  expression  in  Eq.  (8).  Note  that  the  five  coefficients  bs  which  depend  on  the 
actual  A p  are  easily  calculated  because  they  are  just  five  scalars  that  depend  only 
on  A p  and  dp  irrespective  of  the  size  of  the  finite  element  mesh. 

Finally,  the  mass  and  stiffness  matrices  for  the  ith  component  of  the  NX- 
PROMs  can  be  partitioned  as  follows 


M(p0  +  A p)?x 
K(p0  +  A  p)fY 


MA 

MAp,i 

KA„,i 

kcn 

p,i 

KS 

kn 

(11) 


where  superscript  C  indicates  a  constraint  partition,  and  superscript  N  indicates  a 
nominal  mode  partition. 


4.  Local-Interface  Reduction 

If  the  finite  element  mesh  used  is  very  fine,  the  size  of  the  reduced  system-level 
matrices  is  dominated  by  the  constraint  DOFs  corresponding  to  the  C  partition  in 
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Eq.  (11).  The  constraint  DOFs  of  matrices  constructed  by  CMS-based  techniques 
are  difficult  to  reduce.  This  is  an  important  issue  because,  if  the  constraint  DOFs 
cannot  be  reduced,  then  the  overall  structure  cannot  be  efficiently  divided  into 
many  substructures.  To  address  this  issue,  Castanier  et  al.  [6]  suggested  the  use 
of  characteristic  constraint  (CC)  modes.  This  technique  is  based  on  performing 
a  secondary  eigenanalysis  of  the  constraint  partition  (C)  of  the  system-level  mass 
and  stiffness  matrices  constructed  by  CB-CMS.  This  technique  is  applied  after  the 
system-level  matrices  are  constructed.  However,  the  core  idea  of  NX-PROMs  is 
that  all  analyses  are  accomplished  at  the  substructure-level,  and  not  at  the  system- 
level.  Thus,  an  alternate  interface  reduction  technique  is  proposed  next.  The  new 
approach  is  applied  at  the  substructure-level,  and  it  is  referred  to  as  local-interface 
reduction  (FIR). 

The  local-interface  reduction  technique  is  also  based  on  a  secondary  eige¬ 
nanalysis  of  the  constraint  partition.  However,  the  secondary  eigenanalysis  is 
executed  on  the  constraint  partitions  ( C )  of  the  substructural  matrices,  not  the 
system-level  matrices.  The  secondary  eigenanalyses  on  constraint  DOFs  ( C )  of 
the  mass  and  stiffness  matrices  of  component  i  constructed  by  either  CB-CMS  or 
the  NX-PROM  approach  are  given  by 

ii  -  =  0, 

where  A Apti  is  a  diagonal  matrix  which  contains  the  eigenvalues,  and  are  the 

characteristic  constraint  (CC)  modes  of  the  ith  substructure.  They  are  truncated 
for  the  frequency  range  of  interest  by  using  the  eigenvalues  in  Aa p,%.  The  rows  of 
the  CC  modes  indicate  the  constraint  DOFs  of  the  substructure,  and  the  columns 
represent  the  set  of  truncated  CC  modes.  The  CC  modes  for  each  substructure  are 
used  to  reduce  the  interface  DOFs  for  each  boundary  locally.  Note  that  joining  all 
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CC  modes  for  each  interface  between  different  components  may  lead  to  vectors 
which  are  not  necessarily  linearly  independent.  However,  they  span  the  adequate 
space.  Thus,  the  left  singular  values  of  the  CC  interface  modes  may  have  to  be 
used  for  certain  interfaces.  To  demonstrate  the  LIR  procedure,  a  simple  plate 
model  is  used  in  Section  6.2. 

The  set  of  orthogonal  basis  vectors  used  for  all  interfaces  that  a  component  i 
has  with  other  components  are  grouped  in  a  block  diagonal  matrix  which  contains 
the  entire  interface  component  i  has.  The  number  of  blocks  is  equal  to  the  number 
of  components  that  are  connected  to  component  i.  These  matrices  are  denoted 
by  Uap,i.  Of  course,  if  component  i  connects  to  only  one  other  component,  then 
there  is  only  one  block  in  UaP]1.  Next,  the  mass  and  stiffness  matrices  in  Eq.  (1 1) 
are  projected  using  UaRj  as  follows 


Mg,  =  Ug,Mg,UAp,„  Mgf  =  Ug,M™,  Mgg  =  Mg^,, 


KCC  —  TTT  ~KC  TTa  xrCCN  _  tjt  tsNCC  _  p<-CiV  TT  . 

-*'-Ap,i  uAp,i-*vAp,iuAp,D  *^A p,i  u  Apj^Apji  *^A p,i  ^Ap.i u  Ap,i- 


Thus,  the  final  mass  and  stiffness  matrices  with  reduced  constraint  DOFs  are  given 
for  component  i  by 


ivr LIR  — 

iViA  P,i  — 


Mg,  Mgf 
Mgf  Mg, 


1 y-LlR  _ 

)  -^Ap.i  — 


t (CC  tsCCN 

^A  p,i  ^A  p,i 


tsNCC 

^A  p,i 


kn 

^A  p,i 


where  superscript  LIR  indicates  that  the  matrices  are  constructed  using  local- 
interface  reduction.  An  example  is  provided  in  Section  6.2. 


5.  Assembly 

To  predict  the  system-level  dynamics,  the  mass  and  stiffness  matrices  obtained 
in  Sections  3  and  4  for  each  substructure  have  to  be  assembled.  To  do  that,  ge- 


16 


ometric  compatibility  conditions  must  be  enforced.  In  the  following,  we  discuss 
separately  the  case  where  LIR  is  not  used  and  the  case  where  it  is  used. 

Let  us  consider  the  case  where  the  geometric  compatibility  conditions  are  used 
for  models  without  LIR.  In  this  case,  the  constraint  partitions  (C)  of  component- 
level  matrices  keep  the  meaning  of  the  physical  interface  DOFs  of  matrices  ob¬ 
tained  from  FEMs.  This  means  that  the  geometric  compatibility  conditions  be¬ 
tween  interface  DOFs  (constraint  DOFs)  can  be  applied  directly  to  construct  the 
system-level  matrices.  The  complete  component-level  equations  of  motion  for 
component  i  based  on  CB-CMS  or  the  NX-PROM  approach  can  be  expressed  as 

MfOMqfOM  +  KfOMqfOM  =  FfOM,  (13) 


where  ROM  indicates  component-level  matrices  obtained  using  either  CB-CMS 
or  the  NX-PROM  approach.  The  stiffness  matrices  in  Eq.  (13)  obtained  for  com¬ 
ponents  without  parameter  variability  can  be  expressed  as 

r  -i  CBCMS 


ts-R.OM  _  t {CBCMS 


K?  K?n 
Kfc  Kf 


(14) 


For  component  with  parameter  variability,  the  stiffness  matrices  in  Eq.  (13)  can 
be  expressed  as 


t {ROM  =  tsNX 


K 

K 


C  tsCN 
A  p,i  VA  p,i 

NC  T CN 

A  p,i 


NXPROM 


(15) 


The  formulas  for  the  mass  matrices  are  similar  to  those  for  the  system  matrices  in 
Eqs.  (14)  and  (15)  (and  are  omitted  here  for  the  sake  of  brevity). 
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Next,  the  matrices  in  Eq.  (13)  are  grouped  for  all  i  to  obtain 


M 

K 

F 


Bdiag  Mf 0M  •  •  •  M^om 

Bdiag  Kf 0M  •  •  •  Kf OM 

r  i T 

jpi?OMT  . . .  yromT  , 


(16) 


where  n  is  the  number  of  components,  and  Bdiag [•]  denotes  a  block-diagonal 
matrix. 

The  geometric  compatibility  condition  for  the  ROM  is  expressed  as 


(17) 


where,  q,  and  q;  are  the  generalized  coordinates  for  the  constraint  partitions  (C  for 
CB-CMS  or  NX-PROMs)  that  correspond  to  the  interface  between  components 
i  and  j.  Of  course,  there  is  no  compatibility  condition  to  be  enforced  for  two 
components  which  do  not  have  a  common  interface. 

Equation  (17)  is  used  to  transform  the  matrices  in  Eq.  (16)  similar  to  the 
assembly  process  in  all  finite  element  modeling  methods.  The  final  assembled, 
reduced-order,  system-level  matrices  are  given  by 
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where 


K  CN 

K  N 


KCN  kcn  ...  KCN 

Kf  0  •  •  •  0 

0  Kf  0 

0 

o  o  Kf 


K^c  _  KCArT 


and  KG  is  a  matrix  which  is  obtained  by  the  assembly  of  each  interface.  In  gen¬ 
eral,  Kc  is  a  matrix  which  has  a  smaller  size  than  the  C  partitions  in  K.  The  same 
process  is  applied  to  obtain  Ff  .  Also,  similar  relations  are  obtained  for  the  mass 
matrix  M  (and  are  omitted  here  for  the  sake  of  brevity). 

Finally,  the  compatibility  conditions  for  models  constructed  using  LIR  can  be 
expressed  almost  identically  to  those  for  models  without  LIR.  The  only  difference 
is  that  the  generalized  coordinates  in  Eq.  (17)  represent  amplitudes  of  character¬ 
istic  constraint  modes  or  amplitudes  of  the  basis  vectors  used  to  capture  the  space 
spanned  by  the  characteristic  constraint  modes. 


6.  Numerical  Results 

6.1.  Element-level  nonlinearity  of  brick  and  other  types  of  finite  elements 

Figure  (3)  shows  a  simple  plate  structure  modeled  with  a  shell-type  finite  el¬ 
ements,  where  t  is  the  thickness  of  the  plate  (t  =  0.2  mm).  To  examine  the 
variation  in  the  entries  of  the  finite  element  mass  and  stiffness  matrices  for  this 
substructure,  the  thickness  t  is  varied  by  increments  of  At  =  0.01  mm. 

The  size  of  the  mass  and  stiffness  matrices  of  the  simple  structure  shown  in 
Fig.  3  is  15,  582x15,  582.  As  an  example,  the  32nd  diagonal  entries  of  the  mass 
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and  stiffness  matrices,  M32)32  and  K39.39,  are  shown  in  Fig.  4  as  the  thickness 
varies  (approximated  matrices  are  shown  by  a  dash-dot  line,  and  the  actual  ma¬ 
trices  are  shown  by  a  solid-line  for  various  values  of  A p).  One  may  observe  that 
those  entries  of  the  mass  and  stiffness  matrices  vary  almost  linearly.  To  capture 
these  variations  very  accurately,  a  cubic  interpolation  is  used  as  follows 

K(p0  +  A p)  «  K0  +  KiA p  +  K2A p2  +  K3A p3.  (18) 

The  matrices  K0,  Ki,  K2,  and  K3  are  computed  by  using  the  stiffness  matrices 
evaluated  at  three  parameter  values.  These  values  are  the  reference  value  p0,  the 
first  perturbed  value  Pi  =  Po  +  Sp,  the  second  perturbed  value  p2  =  po  +  2 Sp  and 
the  third  perturbed  value  p3  =  p0  +  3 Sp.  The  procedure  is  standard  and  similar  to 
the  one  used  in  Section  2,  so  its  details  are  omitted  here  for  the  sake  of  brevity.  A 
similar  interpolation  is  used  for  the  mass  matrix.  Note  that,  in  contrast  to  Taylor 
series,  the  cubic  interpolation  does  not  require  the  calculation  of  derivative  terms. 

Next,  to  examine  the  variations  in  the  entries  of  the  mass  and  stiffness  matri¬ 
ces  for  brick-type  elements,  the  same  plate  structure  is  modeled  with  brick-type 
elements.  The  nominal  thickness  of  the  plate  is  (the  same)  0.2  mm  and  is  varied 
by  (the  same)  increments  of  At  =  0.01  mm,  as  shown  in  Fig.  5.  The  thickness  is 
varied  in  a  region  near  the  center  of  the  plate.  The  entries  of  the  mass  and  stiffness 
matrices  near  the  DOFs  where  the  thickness  is  varied  are  affected.  A  sample  DOF 
affected  is  the  645th.  The  645th  diagonal  entry  of  the  mass  matrix  varies  linearly 
(and  is  omitted  for  the  sake  of  brevity).  The  same  entry  of  the  stiffness  matrix, 
however,  does  not  vary  linearly,  as  shown  in  Fig.  6  (left),  where  exact  values  are 
shown  by  a  solid-line.  To  capture  this  nonlinear  variation,  the  cubic  interpola¬ 
tion  function  in  Eq.  (18)  is  used.  The  approximate  values  obtained  are  shown  by 
a  dash-dot  line.  These  results  show  that  Eq.  (18)  is  not  good  enough  to  capture 
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the  highly  nonlinear  variation  of  the  stiffness  matrix.  Therefore,  a  fourth-order 
interpolation  is  used  as  follows 

K(p0  +  A p)  «  K0  +  KiAp  +  K2A p2  +  K3Ap3  +  K4Ap4.  (19) 

This  fourth-order  interpolation  captures  well  the  nonlinear  variation  in  the  entries 
of  the  stiffness  matrix  as  shown  in  Fig.  6  (right). 

Based  on  the  results  in  Fig.  6,  one  may  assume  that  the  errors  obtained  based 
on  Eq.  (19)  are  negligible.  However,  that  is  not  correct,  as  demonstrated  by  the 
forced  response  of  the  plate  calculated  using  exact  and  approximated  matrices. 
Figure  (7)  shows  the  forced  response  at  one  of  the  DOFs  on  the  plate  for  excita¬ 
tion  frequencies  near  the  first  resonance.  The  solid-line  represents  the  response 
computed  by  the  actual  mass  and  stiffness  matrices,  and  the  dash-dot  line  in¬ 
dicates  the  response  computed  by  the  mass  and  stiffness  matrices  parametrized 
using  Eq.  (19)  for  the  case  of  At  =  0.37  mm.  Clearly,  the  forced  response  com¬ 
puted  by  the  parametrized  mass  and  stiffness  matrices  does  not  agree  with  that 
computed  by  the  actual  matrices.  This  means  that  the  errors  in  the  entries  of  the 
stiffness  matrix  from  the  fourth-order  interpolation  are  not  small  enough  to  cap¬ 
ture  accurately  the  dynamic  response  of  the  structure  with  a  modified  thickness. 
Note  that,  the  errors  in  the  dynamic  response  are  induced  by  inaccuracies  in  the 
stiffness  matrix,  not  in  the  mass  matrix.  These  results  demonstrated  that  a  new  pa¬ 
rameterization  technique  focused  on  capturing  element-level  nonlinear  variations 
in  the  stiffness  is  needed  for  brick-type  finite  elements. 

The  parametrized  stiffness  matrix  was  calculated  using  Eq.  (8)  for  A p  = 
0.37  mm.  As  a  sample  of  results,  Fig.  8  shows  that  the  645th  diagonal  entry  of 
the  exact  and  the  approximated  stiffness  matrices  match  extremely  well.  Similar 
matches  are  observed  for  all  entries  of  the  matrices.  Next,  forced  responses  were 


21 


calculated  using  these  matrices.  The  results  are  shown  in  Fig.  9.  The  solid-line 
indicates  the  response  computed  using  the  exact  matrices,  and  the  dash-dot  line 
indicates  the  response  computed  by  using  the  new  parametrized  matrices.  The  re¬ 
sults  closely  match.  Moreover,  the  natural  frequencies  for  the  exact  matrices  and 
the  new  parametrized  matrices,  match  also,  as  shown  in  Table  1. 


6.2.  Example  of  local-interface  reduction 

Consider  a  structure  composed  of  5  substructures  and  3  boundaries  as  shown 
in  Fig.  10.  r,  represents  the  jth  boundary.  First,  using  the  constraint  partitions 
(■ C )  of  the  reduced  mass  and  stiffness  matrices,  the  CC  modes  are  computed. 
These  are  cv  X,  <&^3,  <F^'4  and  <E>^;'5.  ,  has  interface  DOFs  for 

boundary  T i  and  T2,  while  has  interface  DOFs  for  T 2  and  r3.  Substructures 

3,  4  and  5  each  have  one  boundary  V2,  and  F3,  respectively.  The  mathematical 
representation  for  these  partitions  for  each  CC  mode  can  be  expressed  as 
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By  using  each  boundary  partition  of  the  CC  modes,  the  augmented  set  of  CC 
modes  for  each  boundary  is  constructed  as 


Pi  _ 

Ap,aug 
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,r2  _ 
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;  (*si)rs 

(sx/2 ; 

,  (20) 

r3 

Ap,aug 

Equation  (20)  describes  the  augmented  CC  bases  for  boundaries  rl5  T2,  and  T3. 
However,  the  augmented  bases  <F^'p  ail„,  ^Ap,aug’  and  ^Ap.aug  are  not  guaranteed 
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to  be  linearly  independent.  Thus,  they  cannot  be  directly  used  to  reduce  the  con¬ 
straint  DOFs  because  of  lack  of  numerical  stability.  Thus,  an  orthogonal  basis 
for  the  space  spanned  by  the  augmented  CC  basis  is  used.  Specifically,  the  left 
singular  vectors  U^'p,  U^p  and  of  the  three  augmented  CC  bases  3>Ap,aug> 
^ Ap.au g’  and  ^Ap.aug  'n  Eq.  (20)  are  computed  for  each  substructure,  and  the  left 
singular  vector  corresponding  to  singular  values  larger  than  0.01%  of  the  maxi¬ 
mum  singular  value  are  selected  for  each  boundary.  The  rows  of  the  orthogonal 
bases  U^,  U^p  and  U^3p  are  (re)sorted  for  each  substructure  to  match  the  inter¬ 
face  DOFs  for  each  boundary.  The  resorted  matrices  are  denoted  by  U^p,  U^p 
and  U|yp.  The  (re)sorted  matrices  are  grouped  for  each  component  i  to  obtain  ma¬ 
trices  UApp  which  UAp,i  are  used  to  project  the  interface  DOFs  onto  the  secondary 
generalized  coordinates  (CC  domain).  For  example,  substructure  1  includes  the 
T i  and  r2  boundaries.  To  reduce  the  interface  DOFs  of  substructure  1,  the  orthog¬ 
onal  bases  U^p  and  U^p  are  grouped  in  a  matrix  Uap,i  given  by 


Uap.i 


u£,  or' 

0r-  V% 


(21) 


The  orthogonal  basis  for  substructure  2  is  constructed  in  the  same  way,  to  obtain 


UAp,2 


u%  or> 
0rs  u% 


(22) 


For  some  of  the  substructures,  the  grouping  of  U^p  matrices  is  not  necessary. 
For  example,  the  bases  used  for  substructures  3,  4  and  5  are  simply  given  by 


Uap,3  =  u£,  UAp,4  =  and  UAp,5  =  Ug,  (23) 


Next,  the  orthogonal  basis  for  each  substructure  is  constructed  using  Eqs.  (21), 
(22)  and  (23),  and  the  interface  DOFs  (C)  of  the  NX-PROMs  generalized  coor- 
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dinates  are  projected  into  the  secondary  generalized  CC  coordinates  as  shown  in 
Eq.  (12). 

6.3.  L-shaped  plate  modeled  with  brick-type  finite  elements 

To  demonstrate  the  proposed  NX-PROM  and  LIR  methodologies,  an  L-shaped 
structure  modeled  with  brick-type  elements  (shown  in  Fig.  11)  and  containing 
thickness  variations  is  investigated  numerically.  Figure  (11)  shows  the  pristine 
structure  and  the  structure  with  thickness  variations.  The  structure  consists  of 
eight  substructures.  Substructures  7  and  8  have  three  cases  of  thickness  variations, 
as  given  in  Table  2.  The  nominal  thickness  of  the  structure  is  1  mm  and  the 
elemental  thickness  is  0.2  mm.  The  thickness  variations  applied  are  very  large 
compared  to  the  nominal  elemental  thickness  considered.  This  causes  the  entries 
of  the  mass  and  stiffness  matrices  to  vary  nonlinearly.  CB-CMS  is  applied  for  the 
1st  to  the  6th  substructure,  and  the  NX-PROM  approach  is  applied  for  the  7th  and 
the  8th  substructures. 

Figure  (12)  shows  the  system-level  forced  responses  of  the  nominal  struc¬ 
ture  and  the  three  cases  of  thickness  variation.  The  dotted  lines  represent  the 
vibration  response  of  the  nominal  structure.  The  crosses  and  circles  represent  the 
responses  of  the  structure  with  thickness  variations  predicted  using  NX-PROMs 
and  full-order  models,  respectively.  To  reveal  the  enhanced  performance  of  the 
NX-PROMs,  the  vibration  response  predicted  by  MC-PROMs  is  also  shown  by 
the  stars  in  Fig.  12.  For  all  three  cases,  the  forced  response  predicted  by  the  full- 
order  models  and  the  NX-PROMs  show  excellent  matching.  However,  the  forced 
responses  predicted  by  MC-PROMs  are  not  accurate  compared  to  the  full-order 
models.  This  is  due  to  the  fact  that  MC-PROMs  cannot  capture  the  volume  varia¬ 
tion  of  brick- type  finite  elements,  which  leads  to  poor  predictions  of  the  vibration 
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response.  Note  also  that  the  thickness  variations  affect  significantly  the  structure, 
as  shown  by  the  significant  differences  between  the  response  of  the  nominal  struc¬ 
ture  and  the  other  responses. 

The  number  of  DOFs  of  the  full-order  model  and  the  NX-PROMs  are  18,300 
and  3,502,  respectively.  The  system-level  DOFs  of  the  NX-PROMs  include  3,000 
interface  DOFs  and  502  generalized  internal  DOFs.  The  number  of  generalized 
internal  DOFs  is  small.  However,  the  number  of  interface  DOFs  is  large,  and 
should  be  reduced.  Thus,  the  LIR  technique  described  in  Section  4  was  applied. 
Fig.  1 1  shows  the  interfaces  which  are  reduced,  where  Ym  indicates  the  mth  inter¬ 
face.  Fig.  13  shows  the  forced  responses  for  the  nominal  structure  and  the  three 
cases  of  thickness  variations.  The  dotted  lines  represent  the  response  of  the  nom¬ 
inal  structure.  The  crosses  and  circles  indicate  the  responses  of  the  structure  with 
thickness  variations  predicted  using  NX-PROMs  with  full  interface  DOFs  and 
full-order  models,  respectively.  The  squares  represent  the  responses  predicted 
using  NX-PROMs  with  LIR.  These  latter  models  use  only  533  interface  DOFs, 
reduced  by  LIR  from  2,498.  The  response  predicted  by  NX-PROMs  with  LIR 
agree  very  well  with  the  responses  of  the  full-order  model. 

6.4.  Results  for  a  high  mobility  multipurpose  wheeled  vehicle  model 

In  this  section,  NX-PROMs  are  used  to  predict  the  dynamic  response  of  a 
realistic  vehicle  model.  We  consider  the  base  frame  of  a  high  mobility  multipur¬ 
pose  wheeled  vehicle  (HMMWV).  The  finite  element  model  for  the  HMMWV  is 
a  conventional  model  used  to  examine  its  dynamic  response  [10,  11,  12,  13,  16]. 
Figure  14  shows  the  system-level  and  substructure-level  finite  element  models  of 
the  HMMWV  frame.  The  cross-bar  structure  is  composed  of  substructures  Cl 
and  Cr ,  which  are  modeled  with  solid-type  finite  elements,  as  shown  in  Fig.  15. 
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The  marked  region  in  Fig.  15  indicates  the  nodes  which  move  due  to  thickness 
variation.  Table  3  indicates  two  cases  of  thickness  variation  for  CL  and  67?.  The 
thickness  variations  applied  are  much  larger  than  those  used  in  the  L-shaped  exam¬ 
ple.  We  chose  these  large  variations  to  demonstrate  that  the  proposed  methods  are 
very  accurate  even  when  the  thickness  variations  are  very  large.  Such  large  vari¬ 
ations  are  encountered  in  practice  especially  when  components  are  re-designed. 
The  structural  and  the  elemental  thicknesses  of  the  67,  and  67, >  substructures  are 
5  mm  and  2.5  mm,  respectively.  NX-PROMs  were  created  for  the  67.  and  67? 
substructures,  and  CB-CMS  was  applied  to  the  remaining  substructures.  Next, 
forces  and  moments  were  applied  to  the  engine  cradle,  and  the  resulting  forced 
response  were  computed.  Figure  16  shows  the  response  of  the  HMMWV  frame 
for  cases  1  and  2.  The  measured  point  is  shown  in  Fig.  14.  The  dotted  line  shows 
the  forced  response  of  the  nominal  HMMWV  structure.  The  circles  and  crosses 
indicate  the  responses  of  the  re-designed  HMMWV  frame  predicted  by  full-order 
models  and  NX-PROMs,  respectively.  The  stars  show  the  responses  predicted  by 
MC-PROMs.  Results  obtained  using  the  full-order  models  and  the  NX-PROMs 
show  excellent  agreement  for  both  cases  1  and  2,  but  the  results  predicted  by  MC- 
PROMs  do  not  agree  well  with  the  response  obtained  using  the  full-order  models. 
Also,  note  that  the  re-design  has  important  effects  on  the  structure,  as  demon¬ 
strated  by  the  significant  difference  between  the  responses  of  the  nominal  and  the 
re-designed  structures. 

The  full-order  model  of  the  HMMWV  has  123,201  DOFs.  The  NX-PROMs 
have  2,683  DOFs,  of  which  1,473  are  constraint  DOFs  and  1,210  are  fixed-interface 
generalized  DOFs.  To  reduce  the  number  of  constraint  DOFs,  LIR  was  applied. 
Figure  14  shows  the  interfaces  of  the  substructures  in  the  HMMWV  frame.  Ta- 
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ble  4  shows  what  are  the  interface  DOFs  for  each  substructure.  Figure  17  shows 
the  response  of  the  HMMWV  model  predicted  by  NX-PROMs  for  cases  1  and  2 
using  different  levels  of  reduction  of  the  overall  interface  DOFs.  Magnified  plots 
near  the  resonant  frequencies  are  included  also.  Note  that  the  accuracy  of  the  re¬ 
sponse  predicted  by  LIR  depends  on  the  number  of  remaining  interface  DOFs.  In 
these  two  cases,  an  acceptable  accuracy  is  obtained  when  the  remaining  interface 
DOFs  are  not  fewer  than  approximately  1,000. 

7.  Conclusions  and  Discussion 

The  key  contributions  of  this  paper  are  as  follows.  The  proposed  next-generation 
parametric  reduced-order  models  (NX-PROMs)  were  developed  by  using  a  new 
parameterization  technique  to  capture  the  element-level  nonlinearity  due  to  vol¬ 
ume  variations  of  finite  elements  of  brick  or  other  types.  In  addition,  to  establish  a 
mathematically  stable  formulation  for  NX-PROMs,  a  new  transformation  matrix 
was  developed  using  a  novel  interpolation  of  static  constraint  modes.  Finally,  a 
local-interface  reduction  (LIR)  technique  was  proposed  for  further  enhancing  the 
computational  efficiency  of  the  NX-PROMs. 

Novel,  next-generation  parametric  reduced  order  models  (NX-PROMs)  for 
predicting  the  vibration  response  of  complex  structures  have  been  presented.  These 
models  address  two  main  drawbacks  of  MC-PROMs.  The  first  is  that  the  param¬ 
eterization  techniques  used  in  MC-PROMs  cannot  capture  the  thickness  variation 
of  brick  and  other  types  of  finite  elements  due  to  element-level  nonlinearity  of 
the  stiffness  matrix.  The  second  drawback  is  that  the  transformation  matrix  for 
MC-PROMs  is  not  numerically  stable.  Thus,  a  new  parameterization  technique 
was  developed  to  capture  the  nonlinearity  of  the  stiffness  matrix,  and  a  new  trans- 
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formation  matrix  was  proposed  to  make  the  NX-PROMs  more  stable  numerically 
and  more  accurate  compared  to  MC-PROMs. 

To  reduce  the  interface  DOFs,  a  new  method  called  local-interface  reduction 
(LIR)  was  developed.  NX-PROMs  were  developed  for  realistic  substructural  anal¬ 
ysis.  In  such  cases,  the  interface  DOFs  should  be  reduced  before  the  system-level 
matrices  are  constructed.  The  LIR  technique  uses  characteristic  constraint  modes 
computed  for  each  substructure  by  using  the  constraint  partition  of  the  reduced 
mass  and  stiffness  matrices  constructed  by  CB-CMS  or  the  NX-PROM  approach. 
By  using  these  characteristic  constraint  modes,  orthogonal  bases  were  defined  to 
reduce  the  interface  DOFs  of  each  substructure.  That  is  a  key  advantage  of  this 
reduction  technique,  and  it  is  very  useful  for  substructural  analysis. 

Similar  to  MC-PROMs,  the  novel  NX-PROMs  also  provide  smaller  system 
matrices  and  shorter  analysis  and  reanalysis  time  to  predict  the  vibration  response 
of  complex  structures.  This  means  that  NX-PROMs  are  especially  useful  for  the 
repetitive  analyses  needed  in  optimization  problems  where  geometric  changes  are 
applied  in  the  design  cycle  for  structures  modeled  with  brick  and  other  types  of 
finite  elements. 
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Table  1:10  lowest  natural  frequencies  for  exact  and  parametrized  matrices  with  volume  variations 


Mode 

Exact 

Approximated 

1 

4278.24 

4278.24 

2 

9024.39 

9011.52 

3 

9068.45 

9053.21 

4 

14323.12 

14323.11 

5 

24952.70 

24952.34 

6 

25024.01 

25023.94 

7 

27463.11 

27460.10 

8 

27656.15 

27654.70 

9 

36098.49 

36098.51 

10 

41275.58 

41262.95 
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Table  2:  Thickness  variations  in  substructures  7  and  8  of  the  L-shaped  plate 


Substructure 

Thickness,  case  1 

Thickness,  case  2 

Thickness,  case  3 

7 

1.00  mm  — >  1.22  mm 

1.00  mm  — y  1.42  mm 

1.00  mm  — y  1.81  mm 

8 

1.00  mm  — >  1.22  mm 

1.00  mm  — >  1.42  mm 

1.00  mm  — >  1.81  mm 

34 


Table  3:  Thickness  variations  in  substructures  Cl  and  Cr  of  the  HMMWV  frame 


Substructure 

Thickness,  case  1 

Thickness,  case  2 

CL 

5  mm  — y  20  mm 

5  mm  32  mm 

cR 

5  mm  — »  20  mm 

5  mm  — >  32  mm 
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Table  4:  Interfaces  between  substructures  in  the  HMMWV  model 
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Global  Coordinates 


Local  Coordinates 


Figure  1:  Sample  8-node  brick  element  with  global  and  local  coordinates 
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Figure  2:  The  case  where  the  value  of  parameter  p  is  between  po+ISp  and  po+(l  +  l)<5p  for  some 
l  value  between  0  and  3 
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Figure  3:  Simple  plate  structure  modeled  with  shell-type  elements 
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Figure  4:  The  32nd  diagonal  entries  of  the  exact  and  the  parametrized  mass  and  stiffness  matrices 
obtained  by  using  a  classic  cubic  interpolation  for  shell-type  elements 
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Figure  5:  Simple  plate  structure  modeled  with  brick-type  elements 
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Figure  6:  The  645th  diagonal  entry  of  the  exact  and  the  parametrized  stiffness  matrices  obtained 
by  using  a  classic  cubic  interpolation  (left)  and  an  classic  fourth  order  interpolation  (right)  for 
brick-type  elements 
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Figure  7:  Forced  responses  calculated  for  A p  =  0.37  mm  using  the  exact  and  parametrized  mass 
and  stiffness  matrices  obtained  based  on  a  fourth-order  interpolation  for  brick-type  elements 
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Figure  8:  The  645th  diagonal  entry  of  the  exact  and  the  parametrized  stiffness  matrices  obtained 
by  using  the  new  interpolation  for  brick-type  elements 
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Figure  9:  Forced  response  calculated  for  A p  =  0.37  mm  using  the  exact  and  the  parametrized 
mass  and  stiffness  matrices  obtained  based  on  the  new  interpolation  for  brick-type  elements 
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Figure  10:  A  simple  structure  used  to  demonstrate  the  local-interface  reduction 
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Figure  1 1 :  An  L-shaped  plate  shown  without  thickness  variation  (left)  and  with  thickness  variation 
(right);  interfaces  between  components  are  denoted  by  F, 
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Figure  12:  Forced  response  predictions  of  the  L-shaped  plate  for  the  nominal  structure  and  for 
cases  1,  2  and  3  computed  using  full-order  models,  MC-PROMs  and  the  novel  NX-PROMs 
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Figure  13:  Forced  response  predictions  of  the  L-shaped  plate  for  the  nominal  structure  and  for 
cases  1,  2  and  3  computed  using  full-order  models,  NX-PROM,  and  NX-PROM  with  LIR 
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Figure  14:  System-level  and  substructure-level  finite  element  models  of  the  frame  of  a  high  mo¬ 
bility  multipurpose  wheeled  vehicle 
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Figure  15:  Nominal  and  re-designed  cross-bar  composed  of  Cl  and  Cr 
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Figure  16:  Forced  response  predictions  for  the  HMMWV  frame  for  the  nominal  structure  and  for 
cases  1  (top)  and  2  (bottom) 
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Figure  17:  Forced  response  predictions  for  the  HMMWV  frame  obtained  using  NX-PROMs  with 
LIR  for  different  levels  of  reduction;  the  total  number  of  DOFs  obtained  for  each  model  using  LIR 
are  indicated  for  case  1  (left)  and  2  (right) 
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